//==============================================================================
// RandomNumberGenerator.cc
//
// Class provides methods for generating random numbers from various
// distributions (normal, etc.).   
//
//==============================================================================
//.......1.........2.........3.........4.........5.........6.........7.........8
//

// Standard library random number generation functions
// 
// rand()   [0,RAND_MAX]
// srand(unsigned int s)
// #define RAND_MAX <integer constant expression >= 32,767> 
// clock_t clock()  

#include "RandomNumberGenerator.hh"

// KJD This is not right.   Forte is picking up some non-forte time function
// I added this typedef just to force a compile, but who knows if it will work. 
typedef long time_t;
#include <time.h>

/*
//------------------------------------------------------
// * SeedWithClockSmall
//------------------------------------------------------
//
void
RandomNumberGenerator::SeedWithClockSmall()
{
  local_time_type seed;

  time(&seed);
  //cout<<"time="<<seed<<endl;
  // Exaggerate the lower order digits
  double fs = (double)seed/(double)10000.0;
  long   ifs = (int) fs;
  fs = fs - (double)ifs;
  ifs = (long)(fs*10000);
  seed=ifs;
  //cout<<"(RND Seed "<<(unsigned int)seed<<") "<<flush;
  srand((unsigned int) seed);
}
*/

//------------------------------------------------------
// * SeedWithClock
//------------------------------------------------------ 
//
void
RandomNumberGenerator::SeedWithClock()
{
  //local_time_type seed;
  long seed; 

  time(&seed);
  //cout<<"time="<<seed<<endl;
  // Exaggerate the lower order digits
  double fs = (double)seed/(double)10000.0;
  long   ifs = (int) fs;
  fs = fs - (double)ifs;
  ifs = (long)(fs*10000);
  seed=ifs;
  //cout<<"(RND Seed "<<(unsigned int)seed<<") "<<flush;
  srandom((unsigned long) seed);
}

//------------------------------------------------------
// * Uniform
//------------------------------------------------------ 
//
double
RandomNumberGenerator::Uniform(double start,double end)
{
  if ((start >= 0.0) && (end >= start)){
	 double Range = end-start;
	 long rNum = random();
	 double Yout = (double)rNum/(double)LONG_MAX;
	 
	 Yout = start+Yout*Range;
	 
	 if ((Yout < start) || (Yout > end)){
		cout<<"start="<<start<<endl;
		cout<<"end = "<<end<<endl;
		cout<<"Range="<<Range<<endl;
		cout<<"rNum = "<<rNum<<endl;
		cout<<"LONG_MAX="<<LONG_MAX<<endl;
	 }
	 return(Yout);
  }
  return(-1.0);
}


//------------------------------------------------------
// * Gaussian
//------------------------------------------------------ 
// Transforms uniformly distributed random generators into 
// a Gaussian (normal) random generator with specified 
// mean and standard deviation. 
//
//
double 
RandomNumberGenerator::Gaussian(double mean,double StdDev)
{
  double Variance = StdDev*StdDev;
  double PI = 4.0 * atan(1.0);

  double B   = (double)random()/(double)LONG_MAX;
  double Xin = (double)random()/(double)LONG_MAX;

  B = 2.0*PI*B;
  
  double R = sqrt(2.0*(Variance)*log(1.0/(1.0-Xin)));

  double Yout = R * cos(B);

  Yout = Yout+mean;
  
  return(Yout);
}
